#########################################################################
#To make this figure, you must first run "make_full_panel.R" first.######
#########################################################################
rm(list = ls())
setwd("~/Desktop/QJPS Replication")
data <- read.csv("all_countries_ideal_points_panel.csv")

mean <- tapply(data$cscore, data$year, function(x) mean(x, na.rm = T))
sd <- tapply(data$cscore, data$year, function(x) sd(x, na.rm = T))
num <- tapply(data$cscore, data$year, function(x) length(x))

#subset to the set of countries who we have data for prior to 1950
names <- unique(data$cow[data$year == 1950])
data_small <- subset(data, cow %in% names)

mean.50 <- tapply(data_small$cscore, data_small$year, function(x) mean(x, na.rm = T))
sd.50 <- tapply(data_small$cscore, data_small$year, function(x) sd(x, na.rm = T))
num.50 <- tapply(data_small$cscore, data_small$year, function(x) length(x))

#subset to the set of countries who we do not have data for prior to 1950
names <- unique(data$cow[data$year == 1950])
names2 <- as.vector(unique(data$cow[!(data$cow %in% names)]))
data_small2 <- data[data$cow %in% names2,]

mean.n50 <- tapply(data_small2$cscore, data_small2$year, function(x) mean(x, na.rm = T))
sd.n50 <- tapply(data_small2$cscore, data_small2$year, function(x) sd(x, na.rm = T))
num.n50 <- tapply(data_small2$cscore, data_small2$year, function(x) length(x))


plot(names(mean), mean, type = "l", lwd = 2, xlab = "Year", ylab = "Average of Scores", ylim = c(-2.2, 1), main = "All Constitutions", axes = F)
axis(1, at = seq(1780, 2010, 20), cex.axis = .8)
axis(2, at = seq(-2.5, 2, .5), las = 2)

plot(names(sd), sd, type = "l", lwd = 2, xlab = "Year", ylab = "Standard Deviation of Scores", ylim = c(.6, 1.2), main = "All Constitutions", axes = F)
axis(1, at = seq(1780, 2010, 20), cex.axis = .8)
axis(2, at = seq(-2.5, 1.2, .1), las = 2)


#Countries that have been around since 1950.
plot(names(mean.50), mean.50, type = "l", lwd = 2, xlab = "Year", ylab = "Average of Scores", ylim = c(-2.2, 1), xlim = c(1950, 2010), main = "Comparing Existing to New States", axes = F)
lines(names(mean.n50), mean.n50, lty = 2, lwd = 2)
axis(1, at = seq(1940, 2010, 10), cex.axis = .8)
axis(2, at = seq(-2.5, 2, .5), las = 2)
text(1970, .3, "Existing States")
text(1990, -.5, "New States")

plot(names(sd.50), sd.50, type = "l", lwd = 2, xlab = "Year", ylab = "Standard Deviation of Scores", ylim = c(0, 1.3), xlim = c(1950, 2010), main = "Comparing Existing to New States", axes = F)
lines(names(sd.n50), sd.n50, lty = 2, lwd = 1)
axis(1, at = seq(1940, 2010, 10), cex.axis = .8)
axis(2, at = seq(0, 1.4, .2), las = 2)
text(1990, 1, "Existing States")
text(1990, .45, "New States")









